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NOMENCLATURE 


A  complete  listing  of  the  subroutine  names  and  their 
functions,  the  input  NAMELISTS,  and  the  input  and  output  vari¬ 
ables  is  contained  in  Sections  III. A  and  III.B.  Other 
quantities  not  defined  there  are  listed  below: 


Text  Computer  Program 


a 

a* 

A* 

ASTAR 

CD 

CD 

Meaning 

speed  of  sound 
critical  speed  of  sound 
throat  area 

discharge  coefficient,  defined 
in  Eq.  (11-32) 


CDI  'CD3  'CDS 


CD  1 , CD2 , CD3 


d 


D 


g,h 


discharge  coefficient  constants, 
Eq.  (11-33) 

distance  between  the  throat 
locations  at  the  inner  and  outer 
walls.  Fig.  II. 1 

equations  of  the  inner  and  outer 
nozzle  wall  contours  in  the  x-y 
coordinate  system.  Fig.  II. 1 


gi/hl'g2'h2  SI , H 1 , G2 , H2 

G,H 

M  M 

M*  MSTAR 
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p/p0  PP0 

R  R 


dimensionless  quantities  defined 
in  Eq.  (11-25) 

equations  of  the  inner  and  outer 
nozzle  wall  contours  in  the  R-Z 
coordinate  system.  Fig.  II. 1 

Mach  number 

ratio  of  the  speed  at  a  point  to 
the  critical  speed  of  sound 

used  to  denote  physical  order  of 
magnitude 

static-to-stagnation  pressure 
ratio  • 

radial  coordinate  in  cylindrical 
coordinate  system,  Fig.  II. 1 
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e 


R  ,  ,R 
el  '  e 


Computer  Program  Meaning 

average  radius  of  curvature  for 
the  two  bounding  walls  at  the 
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rci , rcs  radii  of  curvature  of  the  inner 

and  outer  nozzle  walls  at  the 
nozzle  throat 
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Re 


24 
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U/ V 
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Ul/V  ,U  ,V  ,  U1 , VI , U2  t  V2 , 
U5,V,  U3.V3 
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X,Y 


y,  'Y 

i  • 


radial  coordinates  of  the  inner 
and  outer  wall  throat  locations, 
Fig.  II. 1 

Reynolds  number  based  on  sonic 
conditions  and  twice  the  throat 
separation  distance,  d 

dimensionless  velocity  com¬ 
ponents  in  the  x-y  coordinate 
system  defined  in  Eqs.  (II-5) 
and  (II-6),  Fig.  II. 1 

transonic  perturbation  velocity 
components  defined  in 
Eqs.  (11-10)  and  (11-11) 

transonic  perturbation  velocity 
components  of  the  first  three 
orders  defined  by  the  expansions 
in  Eqs.  (11-19)  and  (11-20) 

velocity  components  in  the 
cylindrical  R-Z  coordinate 
system.  Fig.  II. 1 

rotated  coordinate  system  non- 
dimensionalized  with  respect  to 
the  throat  separation  distance, 
d,  and  oriented  such  that  the 
y-axis  lies  along  the  minimum 
area  cross-section  and  the 
origin  is  on  the  Z-axis  of  sym¬ 
metry,  Eqs.  (II-3)  and  (II-4) 
and  Fig.  II. 1 

y-  coordinates  of  the  inner  and 
outer  throat  wall  locations. 

Fig.  II. 1 
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Text  Computer  Program 


Meaning 


z  z 

Z  z 


Z* 


ZSTAR 


z 


i 


,Z_ 


ZI  ,Z0 


stretched  axial  coordinate 
defined  in  Eq.  (11-25) 

axial  coordinate  in  the  cylin¬ 
drical  coordinate  system. 

Fig.  II. 1 

displacement  of  the  x-y  origin 
from  the  R-Z  origin  along  the 
Z-axis  of  symmetry.  Fig.  II. 1 

axial  coordinates  of  the  inner 
and  outer  throat  wall  locations. 
Fig.  II. 1 


Greek  Symbols 

8  BETA  inclination  angle  of  the  x-axis 

from  the  Z-axis  of  symmetry, 
positive  counterclockwise. 

Fig.  II. 1 


BETA! 


dimensionless  quantity  defined 
in  Eq.  (11-25) 


Y 


G 


specific  heat  ratio  of  the  gas 


e 


EPS 


expansion  variable  defined  in 
Eq.  (11-18) 


n  ETA 

0  THETA 


parameter  in  the  expansion  vari¬ 
able  definition,  Eq.  (11-18) 

angle  of  inclination  of  the 
velocity  vector  from  the  x-axis, 
positive  counterclockwise 


P 


density 


P* 


critical  density 
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I .  INTRODUCTION 

Annular  supersonic  nozzles  constitute  an  integral  part  of 
a  number  of  devices  of  practical  importance.  These  applica¬ 
tions  include  turbofan  bypass  nozzles,  unconventional  propul¬ 
sion  nozzles  such  as  the  spike,  plug,  and  expansion-deflection 
designs,  as  well  as  the  supersonic-supersonic  ejector  and 
axial-flow  aerodynamic  window.  In  order  to  analyze  the  super¬ 
sonic  flowfield  in  these  nozzles  using  either  the  method-of- 
characteristics  or  a  finite  difference  technique,  an  accurate 
initial  value  line  is  required.  One  natural  place  to  start 
these  calculations  is  in  the  throat  region  of  the  nozzles  using 
an  appropriate  analysis  of  the  transonic  flowfield  which  occurs 
there.  Given  this  starting  line,  the  marching-type  computa¬ 
tions  for  the  steady,  supersonic  portion  of  the  flowfield  can 
then  proceed  in  the  streamwise  direction. 

Several  methods  have  been  utilized  to  analyze  transonic 
flow  in  the  throat  region  of  annular  supersonic  nozzles. 

These  methods  include  inverse  techniques  [l,2]t,  series  expan¬ 
sion  methods  [3-5],  time -dependent  numerical  techniques  [6-8], 
the  method  of  integral  relations  [9],  and  type-dependent 
numerical  relaxation  [10] .  However,  each  of  these  previous 
methods  either  applies  only  to  a  specialized  class  of  annular 

f Numbers  in  brackets  refer  to  entries  in  REFERENCES. 
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nozzles,  such  as  those  with  cylindrical  centerbodies  or 
located  a  large  distance  from  the  axis  of  symmetry,  or  the 
numerical  technique  employed  requires  an  inordinate  amount  of 
computer  time  and/or  memory.  In  addition,  the  inverse  tech¬ 
niques  require  iteration  for  the  direct  problem  of  analyzing 
the  flowfield  in  a  nozzle  of  given  contour.  Only  the  recent 
numerical  methods  of  Cline  [8]  and  Brown,  et  al.  [10]  have 
shown  promise  of  analyzing  nozzle  throat  flows  with  reasonable 
amounts  of  computer  time. 

What  is  desired,  therefore,  is  a  direct  method  which  can 
accurately  and  economically  describe  the  transonic  flowfield 
in  the  throat  region  of  a  large  class  of  annular,  supersonic 
nozzles.  Just  such  a  method  was  developed  by  the  present 
authors  in  [11] .  A  series  expansion  technique  similar  to  the 
one  used  originally  by  Hall  [12]  for  axi symmetric  and  planar 
nozzles  was  utilized  to  find  an  approximate  solution  to  the 
inviscid,  irrotational  governing  equations.  This  solution  may 
be  applied  to  a  variety  of  annular  nozzle  configurations 
including  those  for  which  the  centerbody  and  outer  wall  con¬ 
tours  are  both  circular  arcs  or,  alternately,  those  for  which 
one  boundary  is  straight.  In  addition,  the  main  flow  direc¬ 
tion  may  be  either  parallel  or  inclined  with  respect  to  the 
axis  of  symmetry.  In  the  limit  as  the  centerbody  radius 
approaches  zero  and  the  outer  wall,  respectively,  the  solutions 
for  the  simpler  cases  of  axisymmetric  and  planar  nozzles  are 
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obtained  so  that  these  cases  may  also  be  analyzed.  Since  the 
solution  is  of  the  series  expansion  type,  one  of  its  major 
advantages  is  the  speed  and  reliability  of  its  numerical 
implementation,  making  feasible  parametric  studies  and  itera¬ 
tive  calculations. 

The  purpose  of  the  present  report  is  to  describe  the 
FORTRAN  computer  program  TRANNOZ  which  has  been  developed  to 
implement  the  solution  just  discussed.  After  a  brief  summary 
of  the  theoretical  development  of  the  expansion  solution,  the 
computer  code  is  discussed  in  detail  including  descriptions 
of  its  subroutines  and  functions  and  of  its  input  and  output 
variables.  Input  instructions  are  also  given  together  with  a 
sample  input  file  and  the  resulting  output.  A  listing  of  the 
program  is  included  in  the  Appendix. 


II.  THEORY 


A  sketch  of  the  configuration  to  be  analyzed  is  shown  in 
Fig.  II. 1.  It  consists  of  an  annular  supersonic  nozzle  which, 
in  general,  may  be  inclined  with  respect  to  the  axis  of  sym¬ 
metry.  The  R-Z  coordinate  system  is  the  standard  cylindrical 
coordinate  system,  while  the  x-y  coordinate  system  is  rotated 
in  such  a  planner  that  the  y-axis  lies  along  the  cross  section 
of  minimum  area  in  the  nozzle  throat.  The  x-axis  is  perpen¬ 
dicular  to  the  y-axis,  and  the  origin  of  this  coordinate  system 
is  located  on  the  Z-axis  of  symmetry  a  distance  Z*  from  the 
R-Z  origin.  The  angle  6  is  the  inclination  angle  between  the 
x-axis  and  the  Z-axis,  and  d  is  the  distance  in  the  R-Z 
coordinate  system  between  the  inner  and  outer  throat  wall 
locations.  The  coordinates  of  these  last  two  points  as  well 
as  those  of  the  x-y  origin  and  the  equations  of  the  inner  and 
outer  wall  contours  in  the  meridional  plane  are  also  given  in 
both  the  R-Z  and  x-y  coordinate  systems  in  the  figure.  It  is 
to  be  noted  that  for  the  general  case  of  an  inclined,  annular 
nozzle  the  minimum  area  cross  section  does  not  correspond  to 
the  cross  section  of  minimum  distance  between  the  nozzle  walls. 
Because  of  the  radial  factor  involved  in  calculating  the 
annular  area,  the  minimum  area  cross  section  is  located  nearer 
the  axis  of  symmetry  than  the  minimum  distance  cross  section. 
More  will  be  said  about  this  when  subroutine  ARMIN  is  discussed 
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Under  the  assumptions  of  steady ,  inviscid,  ir rotational, 
adiabatic  flow  of  a  perfect  gas,  the  governing  equations  in 
the  cylindrical  coordinate  system  can  be  taken  as  the 
irrotationality  condition  and  the  "gas  dynamic  equation"  [13], 


U  -  v  ■  0 
a  z 


(u2-a2)u.  +  (v2-a2)v_  +  2UVU_  -  ”  0 


(II-l) 


(II-2) 


a  =  speed  of  sound  , 

where  the  subscripts  are  used  to  denote  partial  differentiation 
with  respect  to  Z  and  R.  Transforming  from  the  R-Z  to  the 
rotated  x-y  coordinate  system  where  lengths  are  non-diraension- 
alized  with  respect  to  the  throat  separation  distance,  d,  and 
velocities  with  respect  to  the  critical  speed  of  Bound,  a*. 


x  =  cose  +  5  *in$ 


y  =  -  ■  sing  +  §  cos3 


(II-3) 


(II-4) 


U  V 

— y  cose  +  —f  Sing 

ci  « 


(II-5) 


U  V 

v  *  -  r*  sine  +  — r  cose  , 
a*  a* 

and  using  the  adiabatic  relation. 


£)*  -  - 


the  governing  equations  become 


u  -  v 
y  * 


(II-6) 


(II-7) 


(II-8) 
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The  next  step  in  the  analytical  development  involves  introduc¬ 
tion  of  the  transonic  perturbation  velocity  components,  u  and 
3,  by  the  relations 

u  -  1  +  u  (11-10) 


v  -  v  ,  (11-11) 

which,  when  substituted  into  Eqs.  (II-8)  and  (II.-9)  result  in 

U  —  V  *  0  (11-12) 

y  x 


(- 2a-*2-  Sr  '  r*r  a+a,aa,  +  (Sr  -*’-2  Sr. 4  *  Sr  a’K 


2  Sr  a 


yTT  ***  ”  [v  cosB  +  (l+u)sin  8] 

y  cosB  +  x  sinS  1 

(11-13) 


The  boundary  conditions  for  this  inviscid  analysis  are 
that  the  nozzle  walls  must  be  streamlines.  Taking  y=g(x)  and 
y=h(x)  as  the  equations  for  the  inner  and  outer  wall  contours, 
respectively,  the  boundary  conditions  may  be  written  as 

v(x,g  (x) )  »  (l+u(x,g(x) ) Jg' (x)  (11-14) 

v(x,h(x) )  -  tl+u(x,h(x))]h' (x)  ,  (11-15) 

where  the  prime  is  used  to  denote  differentiation  with  respect 


to  x. 


To  this  point  in  the  development  no  approximations  to 
either  the  governing  partial  differential  equations  or  the 
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boundary  conditions  have  baen  made.  In  order  to  proceed,  an 
expansion  parameter  must  be  defined  so  that  the  perturbation 
velocity  components  can  be  expanded  in  appropriate  series  and 
substituted  into  the  equations  and  boundary  conditions. 

Based  on  the  experience  of  Kliegel  and  Levine  [14]  and 
Thompson  and  Flack  [5] ,  the  expansion  parameter  used  in  this 


investigation  is 


(E0+n)‘ 


(11-16) 


where  R  is  an  average  dimensionless  radius  of  curvature  for 
the  two  bounding  walls.  The  parameter  n  is  included  in  order 
to  improve  the  convergence  properties  of  the  series  solution 
for  nozzles  with  a  small  wall  radius  of  curvature.  For  n>l, 
e  is  less  than  unity  regardless  of  how  small  Rc  may  be. 
Defining  Rg  in  terms  of  the  second  derivatives  of  the  equa¬ 


tions  for  the  wall  contours. 


R  5 
e 


the  definition  of  e  becomes. 


h" (0)  -  g" (0) 
2  +  n  {h“(’BJ"-"girT 


(11-17) 


(11-18) 


This  is  the  definition  actually  used  in  the  evaluation  of  the 
series  solution. 

The  solution  technique  then  proceeds  by  investigating  the 
orders  of  magnitude  of  the  various  terms  in  the  equations  and 
boundary  conditions  and  by  defining  appropriate  0(1)  quantities. 
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Expanding  the  perturbation  velocity  components  u  and  v  as, 
&(z,y)  =  Uj(z,y)e  +  u2(z,y)e2  +  us(z,y)eS  +  ...  (11-19) 

v (z,y)  »  j^H  e]  ^(z^yjc  +  v2(z,y)e2  +  vs(z,y)e*  +  . .  .J  , 

(11-20) 

substituting  into  governing  Eqs.  (11-12)  and  (11-13)  and 
boundary  conditions  (11-14)  and  (11-15) ,  and  gathering  coeffi¬ 
cients  of  like  powers  of  e  results  in  the  formulations  for  the 
various  solution  orders  in  the  expansion  technique.  For  the 
first  order  solution,  the  governing  equations  are 


o 

II 

f|s 

i 

jfl* 

(11-21) 

3u  3v  6, 

_  _ i  x  _ _  4.  J 

+  v» 

.  s  n 

(11-22) 

2ui  3z  3y 

y 

*  V 

with  corresponding  boundary  conditions 

v.  (*'*.)  -  «.  +  9,* 

(11-23) 

M2,y")  *  h. +  v 

« 

(11-24) 

where , 

9,  * 

hi 

5  [l$i]  1/Jc'S/2h'  (0) 

(11-25) 

2g"(0) 

h 

2h" (0) 

g2  ‘  h" (0)  -  gH(0) 

ft2 

"  h"(0)  -  g"  (O')  * 

In  the  definitions  in  (XI 

-25) 

note,  in  particular, 

the  trans- 

formation  from  the  coordinate  x  to  the  stretched,  axial 
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coordinate  i.  The  formulations  for  the  higher  order  solutions 
are  similar  although  they  contain  many  more  terms  than  the 
first  order  formulation. 

Equations  (11-21) -  (11-24)  and  the  corresponding  equations 
for  the  higher  order  solutions  are  the  ones  which  must  be 
solved  in  order  to  obtain  the  (ux  ,vx ) ,  (ua ,va ) ,  (us  ,vs ) ,  ... 
perturbation  velocity  components.  As  discussed  in  [11] ,  the 
series  solution  has  been  carried  to  the  third  order  using  the 
method  of  Hall  [12].  The  result  is  an  approximate,  analytical 
solution  for  the  perturbation  velocities  consisting  of  a 
rather  large  and  algebraically  complicated  set  of  constants 
and  functions.  However,  because  of  the  closed-form  nature  of 
the  solution,  these  quantities  can  be  rapidly  evaluated  in  a 
straightforward  manner  using  a  digital  computer . 

With  the  resulting  expressions  for  the  (ux ,vx ) ,  (ua,va), 
and  (us,v3)  components  determined,  other  quantities  of  interest 
may  also  be  found.  The  series  expansions  for  the  following 
flowfield  variables  are  given  below:  the  velocity  components 
u  and  v  in  the  x-y  coordinate  system;  M*,  the  ratio  of  the 
local  speed  to  the  critical  speed  of  sound;  6,  the  angle  of 
inclination  of  the  velocity  vector  from  the  x-axis;  the  Mach 
number,  M;  and  the  local  static-to-stagnation  pressure  ratio, 
P/P0 t 

u(x,y)  ■l+a»l+uIe+  uaea  +  ute*  +  ...  (11-26) 


1 

2 


I 


I 

I 

i 

«* 


■:k 

•» 


T 

mi 


0 


i 


B 
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v(z,y)  ■  v  »  e~J  [vj  e  +  vae2  +  v5e5  +  . .  7J 

M*(z/y)  «  (u2+v2)  »  1  +  Uje  +  uae2  +  |us  +  vj 


0(z,y)  «  tan"1  (v/u)  «  ej  /  |va  c  +  (Vj-UjVjje2 

+  (vs‘uivi'ujvi+uJ''1  jc3  ♦  ...] 


(11-27) 


e  +  . . . 


(11-28) 


(11-29) 


M(z,y) 


-  1  +  (*£)  U,£  +  f»,  +  I  (Y-llujjc3 
*  [Us  +  <  *  7  (Y-l)u,«a  ♦  ujlt3  +  ... 


1  -  »* 


u*  £s 

8  ui  E 


Uj  e  +  . . . 


(11-30) 


D  r  v-l  2-|Y/(y-0  f  y  )Y/(Y*1)[  r 

£-  (z,y)  -  [l  -  IZj-  M*  |  »  1  -  y[ui£  +  u 

+  (us  +  ^  VJ  “  ^  u*)e*  +  ...]  •  (11-31) 


Another  quantity  of  interest  is  the  discharge  or  flow 
coefficient,  CD,  which  is  defined  as  the  ratio  of  the  actual 
nozzle  mass  flowrate  to  that  obtained  from  the  ideal  approxima¬ 
tion  of  uniform,  sonic  flow  at  the  throat. 


c°-C[^u$],.0  • 


(11-32) 


Substituting  the  appropriate  expressions  for  the  quantities  in 
the  integrand  and  carrying  out  the  integration,  the  relation 


for  CD  becomes 

[Cdi  +  coae  +  cos£a  +  “•]  (11-33) 

where  CDl ,  CDa,  and  CD$  are  constants. 

Thus  with  the  (u^Vj),  (ua,va),  (us,v5)  transonic  pertur 
bation  velocity  components  and  the  CDI ,  CDa ,  CQ3  discharge 
coefficient  constants  determined,  all  of  the  flowfield  vari¬ 
ables  of  interest  are  known  to  third  order  in  the  present 
series  approximations.  For  further  details  concerning  the 
development  of  the  expansion  solution  including  the  extensive 
series  of  checks  and  parametric  studies  which  have  been  per¬ 
formed,  reference  [11]  should  be  consulted. 
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III.  COMPUTER  PROGRAM  TRANNOZ 

The  TRANNOZ  computer  code  is  a  FORTRAN  program  for 
analyzing  the  transonic  flowfield  in  the  throat  region  of 
annular  supersonic  nozzles.  This  task  is  accomplished  by 
evaluating  the  series  expansion  solution  developed  in  [11] 
and  outlined  in  the  preceding  chapter.  Since  the  series  eval¬ 
uation  consists  essentially  of  the  straightforward  calculation 
of  a  set  of  constants  and  functions,  and  since  the  solution 
can  be  applied  to  a  variety  of  nozzle  configurations  including 
the  annular,  planar,  and  axisymmetric  cases,  many  nozzle 
throat  flowfields  of  interest  can  be  calculated  in  a  routine 
and  inexpensive  manner.  This  feature  makes  possible  para¬ 
metric  studies  and  iterative  calculations  as  might  be 
necessary,  for  example,  in  a  design  situation. 

A.  PROGRAM,  FUNCTIONS,  AND  SUBROUTINES 

The  TRANNOZ  code  consists  of  a  main  program,  two  function 
subprograms,  and  twelve  subroutines.  A  brief  description  of 
each  of  these  routines  is  given  below.  A  complete  listing  of 
TRANNOZ  is  contained  in  the  Appendix. 

1 .  PROGRAM  MAIN 

The  main  program  reads  and  writes  the  input  variables, 
calls  subroutines  ARMIN  and  DISCO  to  calculate  some  initial 
parameters  neceetsary  in  the  series  evaluation,  calls  the  four 
main  worker  subroutines  CONTOUR,  STLINE,  XPLANE,  and  ZPLANE  as 
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desired,  and  writes  the  results.  It  Is  to  be  emphasized  that 
under  normal  circumstances  all  reading  and  writing  Is  done  by 
program  MAIM.  If  an  error  condition  Is  encountered,  a 
limited  number  of  diagnostics  are  written  by  subroutine  ERROR 
(see  below) . 

2.  FUNCTION  IBND 

Function  IBND  is  the  equation  of  the  inner  wall  contour 
in  the  explicit  form  R*IBND(Z) .  In  the  present  form  of  IBND, 
the  inner  boundary  may  be  either  a  circular  arc  or  a  straight 
line  in  the  meridional  plane,  although  any  contour  which 
satisfies  the  assumptions  of  the  analysis  could  be  used. 

3.  FUNCTION  OBND 

Function  OBND  is  the  equation  of  the  outer  boundary  in 
the  explicit  form  R-OBND (Z) .  As  for  the  inner  boundary, 
TRANNOZ  allows  the  outer  wall  contour  to  be  either  a  circular 
arc  or  a  straight  line  in  the  meridian  plane,  although  other 
contours  could  also  be  utilized. 

4 .  SUBROUTINE  ARMIN 

Subroutine  ARMIN  locates  the  cross  section  of  minimum 
area  in  the  nozzle  throat  and  calculates  some  initial  param¬ 
eters  necessary  in  evaluating  the  series  solution.  As  dis¬ 
cussed  briefly  in  the  previous  chapter,  for  the  general  case 
of  an  inclined,  annular  nozzle  the  minimum  area  cross  section 
does  not  coincide  with  the  cross  section  of  minimum  distance 
between  the  nozzle  walls,  so  that  the  throat  location  is  not 


known  a  priori.  ARMIN  employs  an  iterative,  numerical  tech¬ 
nique  [5]  to  locate  the  throat  whereby  the  location  of  the 
outer  wall  point  is  first  fixed  and  the  location  of  the  inner 
wall  point  is  varied  until  a  minimum  in  the  cross-sectional 
area  is  found.  The  inner  wall  point  is  then  fixed  and  the 
outer  boundary  point  location  is  varied  until  the  area  is 
again  minimized.  This  process  is  continued  until  the 
fractional  change  in  these  successively  determined  minimum 
areas  is  less  than  10~10.  This  technique  has  been  thoroughly 
tested  by  substituting  the  throat  wall  locations  found 
numerically  into  the  algebraic  equations  resulting  from  the 
constrained  minimization  problemt  of  locating  the  throat  of 
an  annular  nozzle  with  circular  arc  contours.  In  all  cases, 
the  equations  were  found  to  be  satisfied  to  within  a  high 
degree  of  numerical  accuracy  thus  verifying  the  method. 

With  the  throat  location  determined,  ARMIN  then  evaluates 
the  following  geometrical  parameters  shown  in  Fig.  II. 1:  A*, 

d,  Z*,  0,  yi  ,  yo  ,  g'(0),  h*  (0) ,  g"  (0)  ,  and  h"(0)  where  A*  is 
the  throat  area  and  y=g(x)  and  y=h(x)  are  the  equations  of  the 
inner  and  outer  contours  in  the  x-y  coordinate  system.  The 
dimensionless  parameters,  e,  gt  ,  hj ,  g2 ,  h2 ,  and  0j ,  defined 
in  Eqs.  (11-18)  and  (11-25),  are  also  calculated. 

f After  elimination  of  the  two  Lagrange  multipliers,  the  formu¬ 
lation  of  the  constrained  minimization  problem  consists  of 
four  nonlinear,  simultaneous,  algebraic  equations  for  the 
coordinates  of  the  inner  and  outer  throat  wall  points 
(Z  ,R  ),  (Z  ,R  )  in  terms  of  the  coordinates  of  the  centers 
ana  radii  o$  curvature  of  the  circular  arc  boundaries. 
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5 .  SUBROUTINE  DISCO 

Subroutine  DISCO  evaluates  the  discharge  coefficient,  CD, 
for  the  nozzle  under  consideration  by  calling  subroutine 
AATRANS  to  calculate  the  constants  CDl ,  CDJ ,  and  CD3  in  the 
series  approximation  for  CD,  Eg.  (11—33) . 

The  next  four  subroutines  are  called  at  the  user's  dis¬ 
cretion  to  perform  the  major  functions  of  the  program. 

6.  SUBROUTINE  CONTOUR 

Subroutine  CONTOUR  finds  the  R-Z  coordinates  of  the 
points  on  contours  of  constant  Mach  number,  M,  dimensionless 
speed  ratio,  M*,  or  static-to- stagnation  pressure  ratio,  p/p0 • 
A  maximum  of  53  points  is  allowed  on  each  contour.  A  contour 
point  is  found  on  both  the  inner  and  outer  boundaries  and  the 
remaining  points  are  equally  spaced  in  the  y-coordinate  from 
yt  to  yo  .  If  the  contour  points  at  y=yt  and  y=yo  are 
essentially  coincident  with  the  adjacent  points  on  the 
boundaries  (difference  in  y  coordinates  less  than  10-3),  they 
are  omitted.  The  order  of  the  points  on  output  is  from  the 
outer  wall  to  the  inner  and  all  of  the  flow  quantities  along 
the  contour  are  printed. 

7 .  SUBROUTINE  STLINE 

Subroutine  STLINE  calculates  the  flowfield  variables 
along  a  supersonic  initial  value  line  for  starting  method-of- 
characteristics  or  finite  difference  computations  for  the 
supersonic  portion  of  the  flowfield  downstream  from  the  nozzle 


a 
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throat.  The  starting  line  which  is  calculated  is  the  constant 
Mach  number  contour  from  the  throat  wall  location  with  the 
higher  Mach  number.  This  line  is  employed  since,  under  the 
assumptions  of  the  analysis,  it  is  the  most  accurate  of  the 
alternatives  considered.  It  is  realized  that  this  initial 
value  line  may  not  be  compatible  with  the  particular  algorithm 
used  to  analyze  the  supersonic  flowfield.  However,  routines 
to  evaluate  other  starting  lines  can  easily  be  developed  by 
using  STLINE  as  an  example. 

8 .  SUBROUTINE  XPLANE 

Subroutine  XPLANE  evaluates  the  flowfield  variables  of 
interest  at  a  specified  number  of  points  along  planes  of  con¬ 
stant  x-coordinate  (Fig.  II. 1)  in  the  nozzle  throat.  The 
points  are  equally  spaced  in  the  y-coordinate  from  y{  to  yo , 
and  a  maximum  of  51  points  on  each  plane  is  allowed.  This 
subroutine  has  been  used  primarily  to  test  the  series  solution 
obtained  in  [11]  for  the  general  annular  configuration  against 
previous  solutions  for  simpler  geometries. 

9.  SUBROUTINE  ZPLANE 

Subroutine  ZPLANE  computes  the  flowfield  quantities  at  a 
specified  number  of  points  along  planes  of  constant  axial 
coordinate,  Z.  A  maximum  of  51  points  is  allowed  on  each 
plane,  and  they  are  equally  spaced  in  the  radial  coordinate,  R, 
from  the  inner  to  the  outer  wall. 


20 

The  next  six  subroutines  are  secondary  routines  and  are 
called  a  number  of  times  by  the  main  worker  subroutines  in 
order  to  carry  out  specific  tasks. 

10.  SUBROUTINE  TRRZXY 

Subroutine  TRRZXY  transforms  the  coordinates  of  a  point 
from  the  R-Z  cylindrical  coordinate  system  to  the  dimension¬ 
less  x-y  system.  This  transformation  is  the  one  indicated  by 
Eqs.  (II-3)  and  (II-4) . 

11.  SUBROUTINE  TRXYRZ 

Subroutine  TRXYRZ  carries  out  the  coordinate  transforma¬ 
tion  of  a  point  from  the  rotated  x-y  coordinate  system  to  the 
R-Z  cylindrical  system.  This  is  the  inverse  of  the  transfor¬ 
mation  expressed  in  Eqs.  (II-3)  and  (II-4) . 

12.  SUBROUTINE  ITER 

Subroutine  ITER  is  a  general  iteration  subroutine  used  to 
find  the  value  of  the  independent  variable  corresponding  to  a 
given  value  of  the  dependent  variable  in  a  functional  relation¬ 
ship.  The  iterations  are  continued  until  error  tests  on  either 
the  independent  or  dependent  variable  are  satisfied.  This 
routine  is  used  by  subroutine  CONTOUR  for  finding  points  along 
contours  of  constant  M,  M*,  or  p/p0 . 

13 .  SUBROUTINE  VARSOR 

Subroutine  VARSOR  is  utilized  to  determine  which  of  the 
three  dependent  variables  M,  M*,  or  p/p0  is  being  held  constant 
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along  a  contour.  It  is  used  in  conjunction  with  the  iteration 
subroutine  ITER  by  worker  routine  CONTOUR. 

14 .  SUBROUTINE  ERROR 

Subroutine  ERROR  writes  a  limited  number  of  diagnostics 
for  error  conditions  encountered  in  other  subroutines  and 
terminates  program  execution.  Four  of  the  diagnostics  are 
involved  with  the  iterations  in  subroutine  ARMIN  for  locating 
the  minimum  area  section  in  the  nozzle  throat  and  initializing 
the  solution  parameters.  Another  diagnostic  is  written  when 
the  iterations  in  subroutine  CONTOUR  for  a  given  point  on  a 
contour  do  not  converge.  The  final  diagnostic  is  written  when 
subroutine  STLINE  is  called  and  the  Mach  number  at  neither  of 
the  throat  wall  locations  is  supersonic.  Since  none  of  these 
error  conditions  has  ever  been  encountered,  they  will  not  be 
discussed  in  further  detail. 

15.  SUBROUTINE  AATRANS 

Subroutine  AATRANS  is  the  longest  routine  in  the  program. 
It  evaluates  all  of  the  constants,  functions,  and  flowfield 
quantities  in  the  three  term  series  expansion  solution 
developed  in  (11].  In  the  form  presented  here,  all  of  the 
constants  and  functions  in  AATRANS  are  double  precision  vari¬ 
ables  so  that  annular  nozzles  in  the  planar  limit  as  the 
dimensionless  distance  from  the  axis  of  symmetry  becomes  very 
large,  can  also  be  analyzed  without  incurring  significant 
roundoff  errors.  The  penalty  for  adding  this  feature  is 
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increased  compilation  and  execution  times,  but  since  the  pro¬ 
gram  is  extremely  fast  the  effect  of  this  penalty  is  minimal. 


B.  INPUT  AND  OUTPUT  VARIABLES 

As  will  be  discussed  in  more  detail  in  the  next  section, 
input  to  the  TRANNOZ  code  is  achieved  through  the  use  of  six 
NAMELISTS:  PARAM,  CONTROL,  NAMECON,  WAMEST,  NAMEXPL,  and 
NAMEZPL.  A  description  of  each  of  the  input  variables  which 
comprise  these  NAMELISTS  follows. 

1.  NAMELIST  PARAM 

ngeom - an  integer  variable  describing  the  nozzle  geometry. 

Fig.  III.l.  If  NGE0M=1,  both  the  inner  and  outer 
boundaries  are  circular  arcs  in  the  meridional 
plane.  If  ngeom=2,  the  inner  boundary  is  a  straight 
line  and  the  outer  boundary  is  a  circular  arc.  If 
ngeom=3,  the  inner  boundary  is  a  circular  arc  and 
the  outer  boundary  is  a  straight  line.  (default=l) 

A i - Z  coordinate  of  the  center  of  curvature  of  the 

inner  boundary  if  it  is  a  circular  arc  (ngeom=1  or  3) 
or  slope  of  the  inner  boundary  if  it  is  a  straight 
line  (ngeom=2),  Fig.  III.l. 

b i - R  coordinate  of  the  center  of  curvature  of  the  inner 

boundary  if  it  is  a  circular  arc  (ngeom=1  or  3)  or 
intercept  of  the  inner  boundary  if  it  is  a  straight 
line  (ngeom=2) ,  Fig.  III.l. 


rci - radius  of  curvature  of  the  inner  boundary  if  it  is  a 

circular  arc  (NGE0M=1  or  3),  Fig.  III.l. 

as - Z  coordinate  of  the  center  of  curvature  of  the  outer 


boundary  if  it  is  a  circular  arc  (ngeom=1  or  2)  or 
slope  of  the  outer  boundary  if  it  is  a  straight  line 
(NGE0M*3) ,  Fig.  III.l. 

B0 - R  coordinate  of  the  center  of  curvature  of  the  outer 

boundary  if  it  is  a  circular  arc  (ngeom=1  or  2)  or 
intercept  of  the  outer  boundary  if  it  is  a  straight 
line  (NGE0M=3) ,  Fig.  III.l. 
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RC0 - radius  of  curvature  of  the  outer  boundary  if  it  is 

a  circular  arc  (ngeom»1  or  2),  Fig.  III.1 

The  next  four  input  variables  in  NAMELIST  PARAM  are  used 
to  establish  a  "window"  area  in  the  nozzle  throat.  Subroutine 
ARMIN  searches  only  this  region  for  the  minimum  area  cross 
section  and  subroutine  CONTOUR  searches  only  this  region  for 
points  on  the  contours  of  constant  Mach  number,  M*,  or  static- 
to-stagnation  pressure  ratio. 

zimin - minimum  Z  coordinate  on  the  inner  boundary  for 

establishing  the  throat  window.  Fig.  III. 2. 

ZIMAX - maximum  Z  coordinate  on  the  inner  boundary  for 

establishing  the  throat  window.  Fig.  III. 2. 

ZSM IN— -minimum  Z  coordinate  on  the  outer  boundary  for 
establishing  the  throat  window.  Fig.  III. 2. 

ZSMAX-— maximum  Z  coordinate  on  the  outer  boundary  for 
establishing  the  throat  window,  Fig.  III. 2. 

G - specific  heat  ratio  of  the  gas,  y.  (default=1.4) 

eta - parameter  n  in  the  definition  of  the  expansion  vari¬ 

able,  Eqs.  (11-16)  or  (11-18).  (default=2 .0) 

nterm - number  of  terms  from  the  expansion  solution  to  be 

used  in  evaluating  the  nozzle  discharge  coefficient, 
CD.  (default=3) 

2.  NAMELIST  CONTROL 

ncont - an  integer  variable  which  indicates  the  number  of 

contours  of  constant  M,  M*,  or  p/p0  to  be  found, 
(default®^)  NAMELIST  NAMECON  is  read  ncont  times. 

start - a  logical  variable  which  if  .TRUE,  causes  a  super¬ 

sonic  initial  value  line  to  be  found  for  starting 
method-of-characteristics  or  finite  difference  calcu¬ 
lations.  (default®. FALSE.)  If  START®. TRUE., 

NAMELIST  NAMEST  is  read. 
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nxpl - an  integer  variable  which  specifies  the  number  of 

planes  of  constant  x  coordinate  along  which  flow- 
field  quantities  are  to  be  found.  (default*#) 
NAMELIST  NAMEXPL  is  read  nxpl  times. 

nzpl - an  integer  variable  which  indicates  the  number  of 

planes  of  constant  Z  coordinate  along  which  the  flow 
variables  are  to  be  determined.  (default*#) 

NAMELIST  NAMEZPL  is  read  nzpl  times. 

3 .  NAMELIST  NAMECON 

nvar - an  integer  variable  which  determines  which  dependent 

variable  is  to  be  held  constant  along  the  desired 
contour.  For  nvar*1,  the  dependent  variable  is  Mach 
number,  M;  for  nvar=2,  it  is  the  dimensionless  speed 
ratio  M*;  and  for  nvar*3,  it  is  the  static-to- 
stagnation  pressure  ratio,  p/p0 •  (default=l) 

value - the  value  of  the  dependent  variable  along  the 

desired  contour. 

npts - the  number  of  points  to  be  found  along  the  contour 

(minimum=4 ;  maximum*53) . 

nterm— number  of  terms  from  the  series  expansion  solution  to 
be  utilized  in  finding  the  contour. 

4 .  NAMELIST  NAMEST 

npts - number  of  points  to  be  found  along  the  supersonic 

starting  line  (minimum=4j  maximum=53) . 

nterm - number  of  terms  from  the  series  solution  to  be  used 

in  determining  flowfield  quantities  along  the  initial 
value  line. 

5.  NAMELIST  NAMEXPL 

x - x  coordinate  of  the  plane  along  which  the  flowfield 

quantities  are  to  be  determined. 

npts - -number  of  points  along  the  plane  of  constant  x 

coordinate  at  which  the  flow  variables  are  to  be 
found  (minimum*2;  maximum-51) . 

nterm— —number  of  terms  from  the  series  solution  to  be 

utilized  in  finding  the  quantities  of  interest  along 
the  x*constant  plane. 


6 .  NAMELIST  NAMEZPL 


Z — - — Z  coordinate  of  the  plane  along  which  the  flowfield 

quantities  are  to  be  determined. 

NPTS - number  of  points  along  the  plane  of  constant  Z 

coordinate  at  which  the  flow  variables  are  to  be 
found  (minimum9 2;  maximum=51) . 

nterm-— number  of  terms  from  the  expansion  solution  to  be 

used  in  finding  the  quantities  of  interest  along  the 
Z*constant  plane. 


The  first  page  of  output  consists  of  a  listing  of  the 
parameters  from  input  NAMELISTS  PARAM  and  CONTROL ,  as  just 
described,  together  with  the  following  initial  quantities 
determined  by  subroutines  ARMIN  and  DISCO: 


ZI - Z  coordinate  of  the  inner  boundary  throat  location, 

Z.  ,  Fig.  II. 1. 

RI- - R  coordinate  of  the  inner  boundary  throat  location, 

Rj  ,  F  xg .  1 1 . 1 . 

20 - Z  coordinate  of  the  outer  boundary  throat  location, 

Z  ,  F  lg .  1 1 . 1 « 

RS - R  coordinate  of  the  outer  boundary  throat  location, 

Ro  ,  F  xg  .  1 1 . 1 . 

astar - throat  area. 

D— - separation  distance,  d,  between  the  inner  and  outer 

throat  wall  locations  in  the  R-Z  coordinate  system. 
Fig.  II. 1. 

BETA - angle  of  inclination,  6,  of  the  x-axis  from  the 

Z-axis  of  symmetry  (positive  counterclockwise) , 

Fig.  II. 1. 

Y I - y  coordinate  of  the  throat  at  the  inner  boundary,  y 

Fig.  II. 1. 

Y0— - - -y  coordinate  of  the  throat  at  the  outer  boundary,  y 

Fig.  II. 1.  6 


EPS - value  of  the  expansion  variable,  e,  Eq.  (11-18). 

CD - nozzle  discharge  coefficient,  CD,  Eq.  (11-33) . 


The  remaining  pages  of  output  consist  of  listings  of  the 
parameters  in  the  optional  input  NAMELISTS  NAMECON,  NAMEST, 
NAMEXPL,  and  NAMEZPL,  as  used,  as  well  as  the  following  flow- 
field  variables  along  each  contour  or  plane: 


Z— — — axial  coordinate  of  the  contour  point. 

r- - radial  coordinate  of  the  contour  point  or  of  the 

point  on  the  plane  of  constant  Z  coordinate. 

v — — - y  coordinate  of  the  point  on  the  plane  of  constant  x 

coordinate. 

u - - - component  of  velocity,  u,  parallel  to  the  x-axis  non- 

dimensional  i  zed  with  respect  to  the  critical  speed  of 
sound,  Eq.  (II-5)  . 

v - component  of  velocity,  v,  parallel  to  the  y-axis  non- 

dimensionalized  with  respect  to  the  critical  speed  of 
sound,  Eq.  (II-6) . 

M* - dimensionless  ratio  of  the  speed  at  a  point  to  the 

critical  speed  of  sound. 

theta - angle  of  inclination  of  the  velocity  vector  from  the 

x-axis,  0=tan~  (v/u)  (positive  counterclockwise) . 

M - Mach  number,  M,  which  is  the  dimensionless  ratio  of 

the  speed  at  a  point  to  the  speed  of  sound  at  that 
point. 

p/ps— — static-to-stagnation  pressure  ratio  at  a  point. 


C.  INPUT  INSTRUCTIONS  AND  EXAMPLE 

The  input  deck  (file)  is  constructed  in  the  following  man 
ner.  The  first  card  is  a  title  card  which  may  be  any  message 
of  up  to  80  characters  at  the  user's  discretion.  This  message 
can  be  used  for  identification  of  both  the  input  file  and  the 
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output  since  it  is  also  the  first  line  of  output.  The  title 
card  is  followed  by  cards  containing  NAMELISTS  PARAM  and 
CONTROL  where  the  usual  conventions  for  reading  NAMELISTS  are 
observed.  As  discussed  in  the  preceding  section,  the  vari¬ 
ables  in  NAMELIST  PARAM  are  geometrical  and  other  parameters 
necessary  in  the  initialization  of  the  problem,  while  those  in 
NAMELIST  CONTROL  are  variables  which  control  the  further 
operation  of  the  program.  Thus,  these  NAMELISTS  must  always 
appear  in  the  input  deck. 

The  remaining  cards  in  the  input  file  contain,  in  order, 
the  optional  NAMELISTS  NAMECON,  NAMEST,  NAMEXPL,  and  NAME2PL. 
The  number  of  times  each  of  these  NAMELIST  cards  appears  in 
the  input  deck  is  determined  by  the  values  of  the  control 
variables  read  in  NAMELIST  CONTROL.  NAMELIST  NAMECON  is  read 
NCONT  times;  NAMELIST  NAMEST  is  read  if  START=.TRUE. ;  NAMELIST 
NAMEXPL  is  read  nxpl  times;  and  NAMELIST  NAME Z PL  is  read  nzpl 
times.  Note  that  if  any  of  the  control  variables  is  left  at 
its  default  value,  the  corresponding  optional  NAMELIST  does 
not  appear  in  the  input  deck. 

Any  number  of  problems  can  be  solved  with  a  single  input 
file  by  simply  repeating  the  sequence  described  above.  It  is 
important  to  note,  however,  that  the  default  values  are  reset 
at  the  beginning  of  each  new  problem. 

The  dimensional  variables  on  input  are  ai,  bi,  rci,  as. 
80,  RC0,  ZIMIN,  zimax,  Z0MIN,  Z0MAX,  and  z  which  all  have 


dimensions  of  length.  The  units  for  these  input  parameters 
must  be  consistent  and  are  simply  the  units  of  the  R-Z 
coordinate  system.  On  output  these  variables,  together  with 
zi,  ri,  zb,  RS,  and  o,  will  have  these  same  units  and  the 
throat  area,  astar,  will  have  this  unit  squared. 

As  an  example  consider  the  annular  nozzle  shown  in 
Fig.  III. 3.  Both  the  inner  and  outer  boundaries  are  circular 
arcs  in  the  meridional  plane  so  that  ngeom=1.  The  centers  of 
curvature  of  both  surfaces  lie  along  the  Z=0  plane,  ai=a0=0., 
with  radial  coordinates  bi=-1.625  and  BB-2.fi.  The  radii  of 
curvature  of  the  inner  and  outer  boundaries  are  rci=2.0  and 
RC0*1 . J0 ,  respectively. 

The  input  file  for  this  example  is  shown  in  Fig.  III. 4. 
Following  the  title  card,  NAMELIST  PARAM  is  read.  In  addition 
to  the  geometrical  parameters  just  discussed,  the  throat 
"window"  is  set  by  ZIMIN=-/J.5,  zimax=0.5,  Z0MIN=-J0. 5,  and 
Z0MAX«j3.5.  The  remaining  variables  in  PARAM  are  left  at  their 
default  values,  g=1.4,  eta=2.J0,  and  nterm=3.  In  this  example 
five  contours  are  to  be  found  as  well  as  a  supersonic  initial 
value  line.  Thus,  in  NAMELIST  CONTROL  the  values  ncont=5  and 
start*. TRUE,  are  specified.  Variables  nxpl  and  nzpl  remain  at 
their  default  values  of  zero.  The  next  five  cards  are  used  to 
specify  the  values  of  the  parameters  for  each  desired  contour. 
In  each  case  23  points  along  a  constant  Mach  number  contour 
are  desired  with  three  terms  from  the  series  solution  to  be 
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utilized.  Therefore,  nvar=1,  npts=23,  and  nterm*3  are  j 

i 

i 

specified.  The  desired  Mach  number  along  the  first  contour  is  i 

i 

0.6  so  VALUER. 6.  Following  the  NAMELIST  convention,  only  the  \ 

\ 

variables  which  are  changed  need  to  be  specified  on  successive  j 

reads  of  a  given  NAMELIST.  Therefore,  only  the  different  j 

f 

values  of  the  Mach  number,  value*#. 8, 1.0, 1.2,  and  1.4,  are 

i 

included  on  the  next  four  NAMECON  cards.  The  final  card  in 
the  input  file  specifies  that  three  terms  from  the  expansion 
solution  are  to  be  used  to  find  22  points  along  the  starting 
line,  nterm*3  and  npts*22. 

The  corresponding  output  is  shown  in  Fig.  III. 5.  The 
first  page,  as  discussed  in  the  previous  section,  consists  of 
the  title  card,  the  input  variables  from  NAMELISTS  PARAM  and 
CONTROL,  and  some  initialization  parameters  determined  by  sub¬ 
routines  ARMIN  and  DISCO.  The  next  five  pages  contain  listings 
of  the  coordinates  of  the  M=0.6,  0.8,  1.0,  1.2,  and  1.4  con¬ 
tours  as  well  as  the  flowfield  properties  along  these  contours. 

The  variable  nsolv  listed  on  these  pages  is  the  number  of  con¬ 
tour  points  actually  found.  It  is  included  because  the  contour 
may  pass  out  of  the  window  area,  thus  reducing  the  number  of 

points  actually  found,  and  also  because  if  the  contour  points  " 

j 

at  yt  or  yo  are  essentially  coincident  with  those  on  the  | 

corresponding  boundaries,  they  are  eliminated.  The  latter  is  | 

) 

the  case  for  the  M=1.0  and  1.2  contours  for  which  nsolv=22.  » 

The  last  page  of  output  is  a  listing  of  the  supersonic  initial  1 
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value  line,  which  In  this  case  is  the  constant  Mach  number  con¬ 
tour  M*1.16  originating  from  the  outer  throat  wall  location. 

The  Mach  number  contours  listed  in  Fig.  III. 5  are  plotted 
in  Fig.  III. 6  together  with  the  corresponding  data  obtained  in 
[11] .  The  agreement  between  the  data  and  the  series  solution 
is  seen  to  be  quite  satisfactory.  Note  that  the  values  of  R  , 

C  I 

and  listed  in  the  Figure  are  the  radii  of  curvature  of  the 
inner  and  outer  boundaries  non-dimensionalized  with  respect  to 
the  throat  separation  distance,  d,  which  in  this  case  has  a 
value  of  0.625. 

D.  GENERAL  DISCUSSION 

In  addition  to  annular  nozzles,  the  solution  developed  in 
[11]  and  the  TRANNOZ  code  can  be  used  to  analyze  throat  flow- 
fields  for  the  simpler  cases  of  planar  nozzles  and  axi symmetric 
nozzles  without  centerbodies .  However,  some  care  must  be 
exercised  in  specifying  the  nozzle  geometry  for  these  cases  so 
that  the  proper  results  are  obtained. 

The  axisymmetric  nozzle  configuration  is  obtained  in  the 
limit  as  the  inner  boundary  of  the  general  annular  nozzle 
approaches  the  axis  of  symmetry,  yt -*■(),  for  the  case  of  a 
straight  inner  boundary,  NGEOM“2.  The  inner  boundary  cannot 
be  made  to  coincide  with  the  axis  of  symmetry,  yt  *0,  though, 
since  this  leads  both  to  division  by  zero  and  to  zero  as  the 
argument  of  the  natural  log  in  the  evaluation  of  the  various 
constants  and  functions  in  TRANNOZ.  However,  the  y-coordinate 
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of  the  inner  boundary  can  be  made  arbitrarily  small,  e.g., 

mm\  0 

yt =10  ,  thus  providing  the  desired  axisymmetric  results. 

The  planar  configuration,  on  the  other  hand,  is  obtained 
in  the  limit  as  the  dimensionless  distance  from  the  axis  of 
symmetry  to  the  inner  boundary  of  the  annular  nozzle  becomes 
unbounded,  y  -►«>,  since  the  transverse  curvature  effect 
becomes  negligible  in  that  limit.  However,  this  distance  can¬ 
not  be  taken  as  arbitrarily  large  in  running  the  TRANNOZ  code 
because  of  roundoff  error  considerations.  This  is  due  to  the 
fact  that  the  constants  and  functions  in  the  expansion  solu¬ 
tion  are  proportional  to  powers  of  y,  y}  ,  and  yo  so  that  as 
the  latter  quantities  become  large,  the  evaluation  of  the 
desired  quantities  involves  sums  and  differences  of  very  large 
approximately  equal  quantities.  Above  certain  values  of  the  y 
coordinates,  roundoff  error  persists.  Table  III.l  shows  the 
limiting  values  of  yt  above  which  roundoff  error  affects  the 
solutions  for  the  various  orders  of  both  the  perturbation 
velocity  components,  (ut  ,Vj) ,  (u2,v2),  and  (u3»v3)r  and  the  dis 
charge  coefficient  constants,  CD1  ,  Cp2 ,  and  CD3 .  These  limit¬ 
ing  values  are  shown  for  both  single  and  double  precision 
versions  of  subroutine  AATRANS,  the  double  precision  version 
being  the  one  routinely  used  and  the  one  presented  here. 

Notice  that  the  limits  on  the  discharge  coefficient  constants 
are  stricter  than  those  on  the  velocity  components  which  is  a 
result  of  the  definition  of  Cp  as  being  the  integral  of  the 


density-velocity  product,  Eq.  (11-32) .  The  discharge  coeffi¬ 
cient  constants  therefore  contain  higher  powers  of  y,  y.  ,  and 
yo  than  do  the  corresponding  perturbation  velocity  components. 
For  the  double  precision  version  of  AATRANS ,  Table  III.l  shows 
that  for  a  60-bit  machine  and  third  order  solutions,  a  value 
no  larger  than  y  =1000  should  be  utilized  for  investigation  of 
the  velocity  components  and  a  value  no  larger  than  y.  =250 
should  be  used  for  determination  of  the  discharge  coefficient. 
Both  of  these  values  provide  a  very  good  approximation  to  the 
planar  configuration. 

It  is  to  be  noted  that  in  analyzing  the  planar  limit  any 
of  the  three  geometrical  options,  NGEOM=l,  2,  or  3,  can  be 
used,  Fig.  III.l.  In  particular,  utilizing  ngeom=1  allows 
investigation  of  asymmetric,  planar  nozzles,  i.e.,  planar 
nozzles  with  unequal  radii  of  curvature  for  the  two  bounding 
walls. 

Some  limitations  of  the  series  expansion  solution  should 

also  be  mentioned.  During  the  course  of  the  series  solution 

1  /2  3 /2 

development,  the  estimates  x=0(e  '  )  and  tan$/y=0(e  )  were 

made  to  satisfy  order-of-magnitude  consistency  requirements  in 

the  governing  equations.  The  first  estimate  implies  that  the 

solution  is  valid  only  in  the  transonic  region  near  the  throat 

plane,  x=0.  However,  as  demonstrated  in  Fig.  III. 6,  the 

results  appear  to  be  quite  accurate  through  a  wide  region  of 

3  /2 

the  throat.  The  second  estimate,  tan6/y=0(e  ),  means  that 


for  annular  nozzles  located  a  small  dimensionless  distance  from 
the  axis  of  symmetry,  y=0(l),  only  small  angles  of  inclination 
to  the  axis  of  symmetry  can  be  analyzed.  However,  this  would 
seem  to  be  the  only  physically  realistic  case  anyway,  as  one 
would  not  expect  to  encounter  applications  in  which  the  nozzle 
is  both  near  and  highly  inclined  to  the  axis  of  symmetry.  As 
the  distance  from  the  axis  of  symmetry  becomes  larger,  the 
requirement  of  small  inclination  angles  can  be  relaxed.  A 
final  limitation  is  the  result  of  the  definition  of  the  expan¬ 
sion  parameter,  e=(R0*Hi)”1-  As  discussed  more  fully  in  [11], 
the  series  solution  is  limited  to  nozzles  whose  wall  radius  of 
curvature  is  of  the  order  of  the  throat  separation  distance,  d, 
or  larger.  For  nozzles  with  Rc<<d,  unrealistic  results  are 
obtained. 

Despite  these  limitations,  it  is  felt  that  the  TRANNOZ 
code  provides  a  fast,  inexpensive,  and  easy-to-use  tool  for 
analyzing  throat  flowfields  in  a  number  of  nozzle  configurations 
of  interest.  On  a  CDC  7600  computer,  the  compilation  time 
(FTN  compiler)  for  TRANNOZ  is  approximately  11  seconds  while 
the  CPU  time  required  for  the  sample  problem  presented  in 
Figs.  III.4-III.6  is  1.1  seconds.  Also,  since  the  TRANNOZ  code 
was  written  as  a  flexible  subroutine  library,  other  worker  sub¬ 
routines  similar  to  CONTOUR,  STLINE,  XPLANE,  and  ZPLANE  can 
easily  be  developed  to  carry  out  functions  not  currently 
included  in  the  program. 
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Figure  III. 2  Sketch  depicting  input  variables  zimin,  zimax, 
ZSMIN,  Z0MAX 
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EXAMPLE  INPUT  AND  OUTPUT- -MACH  NUMBER  CONTOURS  IN  AN  ANNULAR  NOZZLE 
3PARAM  NGEOM-1 , A1 >0.0,81 «-1 . 625, RCI -2. O, AO-O. 0.B0-2 . 0, RCO- 1 .0, 
ZIMIN--0. 5, ZIMAX-0. 5,Z0MIN--0. 5.Z0MAX-0. 5S 
3 CONTROL  NCONT-5, START® . T . S 
SNAMECON  NVAR* 1 , VALUE® 0. 6, NPTS-23, NTERM-3S 
SNAMECON  VALUE-0. St 
SNAMECON  VALUE- t . OS 
SNAMECON  VALUE- 1.2S 
SNAMECON  VALUE- 1 .4S 
SNAMEST  NPTS-22, NTERM-3S 


Figure  III. 4  Example  input  file 


RESULTS  FROM  THE  TRANNOZ  CODE  FOR  ANALYZING  NOZZLE  THROAT  FLOWS- -BY  J.C.  DUTTON 


Figure  III. 5  Example  output  (cont. 


DATA  FOR  A  CONTOUR  OF  CONSTANT  MACH  NO.  (NVAR*1> 
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Figure  III. 5  Example  output  (cont. 


DATA  FOR  A  CONTOUR  OF  CONSTANT  MACH  NO.  <NVAR=1> 


SUPERSONIC  STARTING  LINE  DATA- -CONSTANT  MACH  NUMBER  CONTOUR 


Table  III.l  Roundoff  Error  Investigation 


(a)  Approximate  values  of  y.  at  which  roundoff  errort  affects  solutions 
for  perturbation  velocity  components: 


Single  Precision 


Double  Precision 


First  Order 


3000 


~10 


Second  Order 

(vv>) 


60 


30,000 


Third  Order 

f  •  V3  1 


15 


1000 


(b)  Approximate  values  of  y  at  which  roundoff  error1  affects  solutions 
for  discharge  coefficient  constants: 


Single  Precision 

Double  Precision 

First  Order 

M 

120 

50,000 

Second  Order 
(C°2  ) 

20 

1500 

Third  Order 
(CD3) 

10 

250 

T0n  the  University  of  Illinois  CDC  Cyber  175  digital  computer 


IV.  CONCLUSIONS 


A  FORTRAN  computer  program,  TRANNOZ,  has  been  developed  to 
analyze  the  transonic  throat  flowfields  in  annular,  planar, 
and  axisymmetr ic  supersonic  nozzles.  The  program  evaluates 
the  series  expansion  solution  developed  by  Dutton  and  Addy  in 
1111 .  Among  its  capabilities  are  options  to  find  contours  of 
constant  Mach  number,  M*,  or  static-to-stagnation  pressure,  to 
calculate  an  accurate  initial  value  line  for  starting  method- 
of -'characteristics  or  finite  difference  calculations,  and  to 
determine  flowfield  quantities  along  various  planes  in  the 
nozzle  throat.  Major  features  of  TRANNOZ  are  its  numerical 
speed  and  reliability  so  that  computations  can  be  carried  out 
easily  and  routinely. 

The  functioning  of  the  various  subroutines  in  the  code 
has  been  described  together  with  the  definitions  of  the  input 
and  output  variables,  detailed  input  instructions,  and  an 
example  input  file  with  the  corresponding  output.  A  brief 
description  of  the  theory  upon  which  the  solution  is  based  has 
also  been  included. 

As  a  result  of  the  characteristics  described  above,  it  is 
felt  that  the  TRANNOZ  program  provides  an  efficient  means  for 
obtaining  approximate  solutions  for  the  flowfields  in  the  throat 
regions  of  a  variety  of  supersonic  nozzles  of  interest. 
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APPENDIX.  TRANNOZ  PROGRAM  LISTING 
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